Collective electromagnetic response of discrete metamaterial systems 
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We develop a general formalism to describe the interactions between a discrete set of plasmonic res- 
onators mediated by the electromagnetic field. The resulting system of equations for closely-spaced 
meta-atoms represents a cooperative metamaterial response and we find that collective interactions 
between asymmetric split-ring resonators arranged in a 2D lattice explains the recent experimen- 
tal observations of system-size dependent narrowing of the transmission resonance linewidth. We 
further show that this cooperative narrowing depends sensitively on the lattice spacing. 
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Resonant multiple scattering plays an important role 
in mesoscopic wave phenomena that can also be reached 
with electromagnetic (EM) fields. In the strong scat- 
tering regime, interference of different scattering paths 
between discrete scatterers can result in, e.g., light lo- 
calization [ij-an effect analogous to the Anderson local- 
ization of electrons in solids. Metamaterials comprise 
artificially structured media of plasmonic resonators in- 
teracting with EM fields. Due to several promising phe- 
nomena, such as the possibility for diffraction-free lenses 
resulting from negative refractive index Q, there has 
been a rapidly increasing interest in fabrication and the- 
oretical modeling of such systems. The role of strong 
scattering and cooperative effects between discrete res- 
onators in metamaterial arrays are largely unknown, as 
such phenomena cannot be easily captured by commonly 
employed theoretical techniques based on the homoge- 
nization approximation of the scattering medium or full 
solutions of Maxwell's equations in small, computation- 
ally accessible, systems. 

In this letter we develop a theoretical formalism of 
collective interactions between plasmonic resonators, or 
meta-atoms, mediated by EM field. In a computation- 
ally efficient model, which also provides physical insight 
into the interaction processes, we treat each meta-atom 
as a discrete scatterer, exhibiting a single mode of cur- 
rent oscillation and possessing appropriate electric and 
magnetic dipole moments. Interactions with the EM 
field then determine the collective interactions within the 
ensemble, resulting in collective resonance frequencies, 
linewidths, and line shifts. The resulting set of equations 
can represent a strong cooperative response in the case 
of sufficiently closely-spaced resonators, necessitating the 
inclusion of interference effects in multiple scattering be- 
tween the meta-atoms. As a specific example, we study 
asymmetric split-ring (ASR) resonators, each consisting 
of two circular arcs of slightly unequal lengths. The cur- 
rents in these ASRs may be excited symmetrically (an- 
tisymmetrically), yielding a net oscillating electric (mag- 
netic) dipole. In a single ASR we derive eigenstates anal- 
ogous to superradiant and subradiant states in a pair 



of atoms and in an ensemble of ASRs we demonstrate 
how strong collective interactions in a discrete metama- 
terial array are responsible for recent experimental obser- 
vations of a dramatic narrowing of the transmission reso- 
nance linewidth (increasing quality factor) with the num- 
ber of resonators Q . Our analysis indicates the necessity 
of accounting for the strong collective response of meta- 
material systems in understanding the dynamics and de- 
sign of novel meta-materials. Strong interactions between 
resonators can find important applications in metamate- 
rial systems, providing, e.g., precise control and manip- 
ulation of EM fields on a sub-wavelength scale Q and 
disorder-related phenomena 

We consider an ensemble of N metamaterial unit ele- 
ments, metamolecules, each formed by n discrete meta- 
atoms, such that the position of the meta-atom j is de- 



noted by Vj (j = 1, 



N). This ensemble is driven 



by an external beam with electric field Ej n (r, t) and mag- 
netic field B; n (r, t) whose frequency components have 
wavelengths much larger than the spatial extent of the 
meta-atoms. We assume that meta-atom j supports a 
single eigenmode of current oscillation whose state can 
be described by a single dynamic variable Qj{t) with 
units of charge and whose spatial profile is described 
by time- independent functions Pj(r) and Wj(r). These 
mode functions are defined such that the polarization 
and magnetization densities associated with atom j are 
Pj(r, t) — Qj(t)pj(r) and Mj(r, t) = Ij(t)\Vj(r) respec- 
tively, where Ij(t) = dQj/dt is the current. We consider 
the Power-Zienau-Woolley Lagrangian in the Coulomb 
gauge for the coupled system of the EM field and mat- 
ter (described by an arbitrary dynamic polarization and 
magnetization densities) [f| 



C = K + C cm + Koui + [Qi$) £ i + T M® 



where C cra is the Lagrangian for the free EM 



(1) 



field, 



V^ ou i is the instantaneous Coulomb interaction, £j(t) = 
J d 3 r E(r, t) ■ pj (r) is an electromotive force on the circuit 
j and 3>j(£) = / d 3 r B(r, t) ■Wj(r) is an effective field flux 
through circuit j, K = Y]- IK /2, and I is a phenomeno- 



2 



logical constant representing an inertial inductance which 
is present in the absence of any EM field interactions. 
This inertial inductance may result, e.g., from the effec- 
tive mass of the charge carriers or surface plasmons as 
they move through the circuit. 

We now derive a set of coupled equations governing the 
excitation of the current oscillations driven by a near res- 
onant field. Following the procedure in Ref. @, we find 
the system Hamiltonian in terms of the dynamic variables 
Qj and their conjugate momenta <pj = 11 j + <frj and arrive 
at the equations of motion Qj(t) = {4>j{t) — $>j(t))/l and 
4>j(t) — £j(t). In order to integrate the scattered fields 
we express the positive frequency components of the 
electric displacement and magnetic induction in terms 
of the normal modes a q , D+(r) = J2q&q£q a q el<ir an d 
B+(r) = V^o/eoEqk x e q ^ q a q e l ^ T , where D(r) = 
D+(r) + D"(r), (D+)* = D~, and q refers to both 
wavevector q and the transverse polarization e q . We ob- 
tain the equations of motion for a q from Eq. |T]). The 
procedure resembles derivation for the expressions of the 
scattered light in an atomic vapor Q and full details will 
be given elsewhere . In the limit that the size of the 
circuit elements is much smaller than the wavelength of 
light, we treat meta-atoms as radiating dipoles and ig- 
nore the higher-order multipolc-ficld interactions. The 
electric and magnetic dipole moment for the j th meta- 



atom are cL 



}jhjdj and m , 



IjAjihj respectively, 



where hj and Aj have units of length and area, defined 
such that hjdj = J d 3 rpj(r) and Ajthj = J d 3 rwj(r). 
The total electric and magnetic fields are the sum of 
the incident field and the scattered field radiated by the 
meta-atoms E(r, t) = E in (r,t) +Es(r, t) and H(r, f) = 
H in (r,f) +H s (r,t), where E s (r,t) = V, E Stj (r,t) and 
Es,j(r, f) denotes the field emitted by the meta-atom j. 
The Fourier components of Eg (r, t) and Hg (r, t) for the 
frequency f2 read (k = fl/c) 



E+(r,0) 



fc 3 
Atteq 



fc 3 



G(r-r',fc)P + (r',ft) 



-G x (r-r',fc)M+(r',ft) 



H+(r,m = — I d 3 r' G(r-r',fc)M(r',f2) 

47T ' 



G(r,fc) = 



cG x (r- r',fc)P+(r',ft) 



2 --(l)/,^ , f rr M uW, 



15(kr) . 



(2) 



(3) 



(4) 



Here P + (r) = P+ (r) , where in the dipole approxima- 
tion Pj"(r) ~ d^~(5(r — rj), and G denotes the radiation 
kernel representing the electric (magnetic) field emitted 
from an electric (magnetic) dipole Q, hn" 1 is the outward 
propagating spherical Hankel function of order n, and 



G x (r, fc) is the radiation kernel for the magnetic (elec- 
tric) field of an electric (magnetic) dipole source. Specif- 
ically, G x acting on a vector v yields (f = r/r) 

■G x (r,fc)-v= — (l-J_)* x v. (5) 

r ikr 

In our model the different circuit elements interact 
via dipole radiation. Each element j also interacts with 
its own self-generated radiation, resulting in an effective 
capacitance Cj and self-inductance Lj associated with 
the circuit current mode. Calculation of Cj and Lj is 
analogous to the derivation of the radiative lincwidth 
of an atom and the values of C, and L,- 



depend on 

the geometric structure of the mode functions pj and 
Wj [8(. Consequently, an isolated meta-atom exhibits 
a resonance frequency ojj = 1/ yjLjCj and the charge 
oscillation results in radiative damping rates for elec- 
tric and magnetic dipoles, TE.j(k) — h'jojjCjk 3 /(67reo) 
and TM.j{k) = ^o^jA^k 3 /(6nLj), respectively. We fur- 
ther simplify by neglecting the frequency dependence of 
V E j M j (the Markov approximation), the inertial induc- 
tance (taking the limit l/Lj — > 0), and assume Tmj <C 



so that a single meta-atom behaves as a simple LC 
circuit with resonance frequency uij and a radiative dissi- 



pation rate Te, 



M,j 



In addition, meta-atoms exhibit 



nonradiative decay To, e.g., due to ohmic losses. 

By including the inter-circuit interactions medi- 
ated by the field we obtain the coupled equations 
of motion for the vector of dynamic variables Q = 
(Qi, Q2, ■ ■ ■ , Qu*n) T and their conjugate momenta <f> = 

(01,02, ■ ■ ■ , 0n*7v) T 



if20=L5 (cj - ir E ) - 3w5^ EW -2/2 C"5Q 



(6a) 



- 3^x0/2 - Sin, 

where we defined diagonal matrices Cjji 



(6b) 



CjSj 



^E,j(k)Sj t j' , and 



The matrix £ = hi [I + iToj- 1 - 
(3/2)w _ 2 g M oj~2]Lz appearing in Eq. (|6a| is the induc- 
tance matrix. The driving terms are given in terms 
of the incident fields as & m ,j = ftjEinfo, fi) • dj and 



AjB t 



(rj,Q) ■ rhj. The radiative interactions 
between the scatterers are described by £7e,m,x- These 
have vanishing diagonal elements and off-diagonal ele- 
ments (j ^ j') given by 



[ft 



M 



[ft 
[^ 



'3,3 



3,3 



(7) 

(8) 
(9) 

The first lines of Eqs. ((Ui,b) are analogous to frequency 
dependent inductance and capacitance matrices, respec- 
tively, but due to the long-range interactions mediated 



'3,3 




3 



by the radiation field, these can substantially differ from 
the quasi-static expressions for which Q x is also absent. 

We consider an incident EM field with a dominant fre- 
quency n and bandwidth Sfl, and introduce the normal 
variables for the meta-atoms whose resonance frequencies 
are centered around ujq 




(10) 



When the EM field and the meta-atoms themselves have 
frequency components in a narrow bandwidth about wo, 
such that 5Cl,T E /Mj, \ujj — ojq\ <^ loq, we may further 
simplify by substituting ujq /c for k in the expressions for 
r E / M and Ge/m/x- Then, from Eqs. ((6|) and (fTOjl . we ar- 
rive at coupled equations for bj's. We make the rotating 
wave approximation, in which we neglect terms oscillat- 
ing at high frequencies ~ 2wo, and obtain the following 
closed system for the current oscillation dynamics 



b(f) =Cb(t) + F(i), 



(11) 



where the coupling matrix and the driving are given by 



iS ~ 2 + I ^ E 



o iu> t 



F =i- 



1 



M 



1 



(12) 
(13) 



Here the detuning 6 = lo — ojq and T = Te + Tu + To , 
where we phcnomcnologically account for ohmic losses 
through To- In the limits we've considered, there ex- 
ist as many distinct collective eigenmodes of oscillation 
as there are meta-atoms in the system. Each collective 
mode corresponds to an eigenvector of the matrix C and 
has a distinct resonance frequency and decay rate, deter- 
mined by the imaginary and real parts of the correspond- 
ing eigenvalue. Eqs. f| 1 1H 13[) resemble the set of equations 
describing a cooperative response of atomic gases to res- 
onant light in which case the coupling is due to electric 
dipole radiation alone @, • 

The crucial component of the strong cooperative re- 
sponse of closely-spaced scatterers are recurrent scatter- 
ing events 0, [H| - in which a wave is scattered more than 
once by the same dipole. Such processes cannot generally 
be modeled by the continuous medium electrodynamics, 
necessitating the meta-atoms to be treated as discrete 
scatterers. An approximate calculation of local field cor- 
rections in a magnetodielectric medium of discrete scat- 
terers was performed in Ref. [12j where the translational 
symmetry of an infinite lattice simplifies the response. 

We next apply the formalism to an array of ASRs. 
A single ASR (metamolecule) consists of two separate 
coplanar circular arcs (meta-atoms), labeled by j <E {I, r} 
and separated by a = r> — r; . The current oscillations in 
each meta-atom produce electric dipoles with orientation 
d (d 1 a) and magnetic dipoles with opposite orienta- 
tions m r = rh; (m r _L a, d). Each element of the ASR 



possesses decay rates Te/m/oj = Tb/az/o- An asymme- 
try between the rings, e.g., resulting from a difference in 
arc length, manifests itself as a difference in resonance 
frequencies with lo t = loq + 5lo and uji = ujq — 8u. The 
interaction matrix C for a single ASR is thus given by 



-iSto - T/2 3if G/4 - 3rS/4 \ 
hfG-3r\S/4 i6u-T/2 ) 



(14) 



where T = y/T E T M , T = T E - T M , G = d • G(a, k) • d, 
and S = d • G x (a, k) ■ m r (which is here assumed to be 
real). To analyze the collective modes of the ASR, we 
consider the dynamics of symmetric c + and antisymmet- 
ric c_ modes of oscillation defined by c± = (b r ± 6;)/v2 
that represent the eigenmodes of the ASR in the absence 
of asymmetry 5oj = 0, From Eqs. (fr2"]) and (fTl"]) . one 
finds c± = (-7±/2 T iA)c± - iSuJc T + (F r ± F t )/y/2, 
where the decay rates j± = T E [1 ± 3Im(G)/2] + Tj\/[1 qp 
3Im(G)/2] + f Re(5) + r , and the frequency shift A = 
-3Rc(G)f/2 - 3rim(5)/4. Excitations of these modes 
possess respective net electric and magnetic dipoles and 
will thus be referred to electric and magnetic dipole ex- 
citations. The split ring asymmetry 8u> introduces 
an EM coupling between these modes: the incident field 
with Ej n || d and Bi„ JL rh only excites the electric mode 
when Suj = 0, but for du ^ it can resonantly pump the 
anti-symmetric magnetic mode. 

As an example of a collective response of a meta- 
material array, we study an ensemble of identical ASR 
metamolecules (with a = ae x and d = de y ) arranged 
in a 2D square lattice within a circle of radius r c and 
with lattice spacing u and lattice vectors Ui = ue x and 
U2 = ue y . The sample is illuminated by a cw plane 



wave E+(r) 



2 £e. y < 



lk r with k = fee. 



coupling to the 
electric dipole moments of the ASRs. A transmission 
resonance through such a sheet was experimentally mea- 
sured in Ref. [3j where the number of active ASRs N 
was controlled by decoupling the ASRs with r > r c from 
the rest of the system with approximately circular shaped 
metal masks with varying radii r c . The resonance quality 
factor increased with the total number of active ASRs, 
saturating at about N = 700. 

The incident field drives all metamolecules uniformly, 
and therefore is phase-matched to collective modes in 
which the electric and/or magnetic dipoles oscillate in 
phase. In the absence of a split-ring asymmetry, only 
modes involving oscillating electric dipoles can be driven. 
These modes strongly emit into the ±e z directions en- 
hancing incident wave reflection. The magnetic dipoles, 
on the other hand, dominantly radiate into EM field 
modes within the ASR plane. This radiation may be- 
come trapped through recurrent scattering processes in 
the array, representing modes with suppressed emission 
rates and reflectance, resulting in a transmission reso- 
nance. Consequently, we study the radiation properties 
of the magnetic eigenmode v m of C which maximizes the 
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FIG. 1: (Color online) The cooperative metamaterial response 
displaying collective resonance narrowing, (a) The resonance 
linewidth 7 of the collective magnetic mode v m in the units 
of an isolated meta-atom linewidth Y (for Te = Vm) as a 
function of the number of metamolecules N; the lattice spac- 
ings u = 1/4A (solid line), 3/8 A (dashed line), 1/2 A (dot- 
dashed line), and 3/2 A (dotted line). The intermediate dot- 
dashed line corresponds to an asymmetry Suj = O.ir, while 
Suj — for all other curves. The lower dot-dashed line incor- 
porates ohmic loss (thermal absorption) To = 0.01X with all 
other curves assuming zero absorption, (b) Shows the over- 
lap O m (bF) of v m with the state b/ excited by an incident 
field resonant on the mode v m (for u = 1/2A, and To = 0). 
For 7 C 1 there is a range of asymmetries Suj for which the 
incident field almost exclusively excites the mode v m . 



overlap O m (b J 4) 



i b A| 2 /I] l |vfb A | 2 with the pure 



magnetic excitation = (1, — 1, • • • , 1, — l) T /\f2N. We 
then show that the introduction of an asymmetry allows 
the excitation of v m by the incident field. This mode 
closely resembles that responsible for the experimentally 
observed transmission resonance [H, EH • 

FigureQJi shows dependence of the resonance linewidth 
7 of v m on the number of metamolecules N for differ- 
ent lattice spacings u. In the absence of ohmic losses 
and for sufficiently small u and Suj, 7 oc 1/N for large 
N. The split ring asymmetry only weakly affects v m . 
For Suj = O.Oir, the curves representing 7 are indistin- 
guishable from those for Suj = 0. For the relatively large 
Suj = OAT, however, 7 is increased for TV > 200. The 
cooperative response is very sensitive to u, resulting in 
different linewidth narrowing. For a larger u = 3/2 A, 
however, 7 is insensitive to N indicating the limit of in- 
dependent scattering of isolated metamolecules and a di- 
minished role of cooperative effects. 

As with an isolated ASR, an asymmetry Suj 7^ gen- 
erates an effective two-photon coupling between collec- 
tive electric and magnetic modes. This is analogous to 
electromagnetically induced transparency: one leg of the 
transition is provided by the coupling of the incident 
field to collective electric dipole excitations, while the 
other is provided by the effective coupling introduced by 
the asymmetry. We illustrate this in Fig. QJb), showing 
the relative population O m (b/) of the collective magnetic 
mode v m produced by a resonant field as a function of 
7 and Suj, where bf is the state induced by the uniform 



driving resonant on v m . One sees that for 7 <C T and 
5lo > 7, one can excite a state in which more than 98% 
of the energy is in the target mode v m . For Suj -C 7, any 
excitation that ends up in v m is radiated away before it 
can accumulate; the array behaves as a collection of sym- 
metric metamolecules. For larger Suj, the population of 
v m decreases since the increased strength of the effective 
two-photon transition begins to excite other modes with 
nearby resonance frequencies. Although the density of 
modes which may be excited increases linearly with N, 
the corresponding reduction of 7 means that a smaller 
Suj is needed to excite the target mode, and there is a 
range of asymmetries for which v m is excited. The nar- 
rowing in 7 combined with the near exclusive excitation 
of this mode implies that for larger arrays the radiation 
from the sheet is suppressed and hence the transmission 
enhanced as in Ref . Q . The observed saturation 0] , how- 
ever, may result from a combination of a fixed Suj, which 
leads to the population of other modes for larger arrays, 
and ohmic losses in the resonators which put a floor on 
the collective decay rate 7. Fig. [T^, shows that an ohmic 
loss rate of T/100 results in the expected saturation of 
quality factor with N. 

In conclusion, we developed a computationally efficient 
formalism describing the collective interactions between 
discrete resonators. In principle 1 , one- can calculate the 
EM response also by having knowledge of the material 
comprising the circuit elements and numerically solving 
Maxwell's equations with a numerical mesh small enough 
to resolve the features of each meta-atom. This, how- 
ever, becomes computationally intractable when the sys- 
tem contains more than a few resonators. In ASRs we 
showed how an asymmetry leads to excitation of collec- 
tive magnetic modes by a field which does not couple 
directly to ASR magnetic moments and results in cooper- 
ative response exhibiting a dramatic resonance linewidth 
narrowing, consistently with experimental findings 
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